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Abstract 

In microfluidics, varying wetting properties, expressed in 
terms of the local slip length, can be used to influence 
the flow of a liquid through a device. We study flow past 
surfaces on which the slip length is modulated in stripes. 
We find that the effective slip length for such a flow can 
be expressed as a function of the individual slip lengths 
on the stripes. The angle dependence of the effective slip 
is in excellent agreement with a recent theory describing 
the slip length as a tensorial quantity. This tensorial nature 
allows to induce a transverse flow, which can be used in 
micro mixers to drive a vortex. In our simulations of a flow 
through a square channel with patterned surfaces we see 
a homogeneous rotation about the direction of the flow. We 
investigate the influence of patterns of cosine shaped varying 
local slip on the flow field depending on the orientation 
of the pattern and find the largest effective slip length for 
periods of stripes parallel to the flow direction. 

1. Introduction 

Microdevices tend to behave very differently from devices 
that are used in everyday life. Surface tension effects are 
dominant at these scales, and micropumps and micro valves 
have been fabricated taking advantage of this principle Q]). 
Inherent in these technologies is the need to develop the 
fundamental science of small devices. 
To understand the fluid dynamics in microflows, the clas- 
sical Navier- Stokes equations cannot be applied due to the 
breakdown of the continuum assumption. Recent computer 
simulations apply molecular dynamics®. However, these 
simulations are usually limited to some tens of thousands 
of particles, length scales of nanometers and time scales 
of nanoseconds ( 3 ). To overcome these shortcomings, meso- 
scopic simulation methods such as the Lattice Boltzmann 
method (LBM) are more applicable and have been rec- 
ognized as a promising approach for simulation of mi- 
croflows ©. It has been shown that in the macroscopic 
limit, the Navier Stokes equations can be recovered© [6]). 



However, most previous LBE models virtually correspond 
to the Navier- Stokes equations on the length scale of the 
lattice constant, and when these models are applied to near- 
continuum flows, the slip boundary conditions must be 
perfectly specified to capture possible boundary effects on 
the flow behavior. Thus for the simulation of fluid flows with 
the LBE method, the development of slip-flow boundary 
conditions has been a critical issue. 

In this paper, we use a recently proposed velocity boundary 
condition for the LBE that is independent of the relaxation 
process during collision and contains no artificial slip © and 
thus simulate slip flows by the use of a parameter ( that 
defines the surface tendency to cause slip. Thus we study 
flow past surfaces on which the slip length is modulated in 
stripes. 

2. Simulation Method 

2.1. Lattice Boltzmann Method 

The lattice Boltzmann method (LBM) is a numerical method 
to solve the Boltzmann equation Eq. ([T]) on a discrete 
lattice (4). The Boltzmann equation expresses how the prob- 
ability /(x, v, t) of finding a particle with velocity v at a 
position x and at time t evolves with time: 

^=vV x / + F-V p /+|J = ^(/), (1) 

where F denotes an external body force, V X , P the gradient in 
position and momentum space, and &(f) denotes the BGK 
collision-operator ©. 

Our simulation takes place on a 3 -dimensional lattice. We 
use a D3Q19 lattice with the vectors Ci pointing to the vari- 
ous neighboring sites according to the index i as explained in 
our previous work (7 ). Note that we express all quantities in 
lattice units, i.e., time is measured in units of update intervals 
and length is measured in units of the lattice constant. 



2.2. Boundary Conditions 



3. Results and Discussion 



On the wall boundary nodes, the distribution function as- 
signed to vectors Ci pointing out of the lattice move out 
of the computational domain in the propagation step, and 
the ones assigned to the opposing vectors are undetermined 
because there are no nodes which the distributions could 
come from. Therefore, on the boundary nodes, special rules 
have to be applied. 

Following ( 9 ), we interpolate by means of linear combination 
between the above two boundary conditions to study slip 
flows(7). The respective weight of the no-slip contribution 
and the full slip contribution is controlled through a slip 
parameter ( by the following general relation: 

/r = a fs + (i - o/r , o 

where ff s denotes the density probability function for flow 
with the partial-slip condition, f™ s denotes the density 
probability function of the no- slip boundary condition, and 
f( s stands for the specular reflection boundary condition 
from which full- slip condition may be obtained. 
For the no-slip boundary condition, we use the formulation 
presented recently by Hecht and Harting (10 ). This boundary 
condition is an explicit local on-site second order flux 
boundary condition for 3D LB simulations on a D3Q19 
lattice, which is independent of the relaxation process during 
collision and contains no artificial slip. For full- slip bound- 
ary condition, specular reflection, as suggested in Ref. ©, 
is used which implies no tangential momentum transfer to 
the wall. 

The slip parameter is defined such that at a value of zero, 
the no-slip boundary condition is applied, while at a value 
of one, the full-slip boundary condition is applied. Thus, 
slip flows can be easily simulated by adjusting the value of 
( between zero and one according to physical parameters 
such as the roughness and hydrophobicity of the surface. 
In our previous work (7 ) we have shown that the partial 
slip boundary condition Eq. ^ shows a slip length which 
is independent of the shear rate. The slip length, measured 
in lattice units and obtained at a given value of the slip 
parameter, is independent of the channel width and of 
the density of the fluid. It depends linearly on the BGK 
relaxation time r of the LB simulation. 
By means of the slip parameter ( the slip length b can be 
adjusted in our simulation. The slip length measured in units 
of the lattice constant a is given by 

b = Tt ( 3 ) 
a 3(1-0' 

which matches the numerical results with a relative error of 
less than 0.03%, or 10 -3 lattice units, on the entire range of 
(. The slip diverges for ( — > 1 because this is the full slip 
case, and it approaches zero in the no- slip case for ( — > 0. 



Simulation Setup:. The simulation volume consists of a 
cubic box of 32 nodes in the x, y, and z-directions. Periodic 
boundary conditions are set up along the y and z -directions. 
We use the BGK relaxation time r = 1 throughout the 
paper. The boundary conditions given in the previous section 
are applied in ^-direction. A constant accelerating force is 
applied to the whole domain. The x-component of the force 
is always zero, but the force can be turned continuously 
in the yz-plane. The setup is relaxed up to 70,000 lattice 
updates until a steady state was reached. The profile of 
fluid flow between two plates is studied in three dimensions, 
where the slip parameter ( is a function of the position on 
the surface, which reflects the pattern as described below. 

3.1. Textured Walls with alternating Stripes 

In our simulations we study the flow past patterned sur- 
faces. First, we restrict ourselves to the case of stripes of 
equal width with alternating slip parameters £i and ( 2 . The 
effective slip length of the flow through such a channel is 
found to be a third order function of the slip length obtained 
from homogeneous walls, given by the form: 

b m = A.(&i + &2)+B.(6? + &|) 

+ C-b 1 -b 2 + D-(bi-b 2 + bl-b 1 ) (4) 
+ E-(bl + bl) 

where b m is the slip length of the modified wall with stripes, 
bi and b 2 are the slip lengths obtained for the homogeneous 
walls of the same slip parameters £i and ( 2 , respectively, 
from Eq. ([3]). A, B, C, D, and E are prefactors, which do not 
depend on the slip parameters. However, they still depend on 
the orientation of the stripes on the wall and on the width of 
the stripes of which the pattern is formed. Note, that Eq. ^ 
is chosen to be symmetric in b\ and b 2 which reflects the 
fact that our stripes all have equal width. 
Throughout our simulations, we use stripes of 4 lattice units 
as width, and we observe a perfect third order fit for Eq. ^ 
as seen in Fig.[T] 

For the fit parameters A, B, C, D, and E in Eq. Q, 
depending on the orientation of the stripes to the direction 
of the accelerating force causing the flow, we find the values 
listed in Tab.Q] 

The flow is driven by an accelerating force which is always 
applied in the whole domain. Consistent with the work by 
Feuillebois et al. (fTTT) . we find that the maximum slip length 
is obtained when the stripes are aligned parallel to this force. 
The dependence of the slip length for walls with stripes 
oriented at the angle d between the direction of the force and 
the extension of the stripes has been previously calculated 




Figure 1. Dependence of the effective slip length on the slip parameter of the stripes, according to the stripe 
orientation being a) parallel and b) perpendicular to the direction of the accelerating force causing the flow. This 
figure shows that Eq. {4} gives a nearly perfect fit to the simulation data obtained as long as the local slip length is 
small compared to the width of the stripes. 



Orientation 


A 


B 


C 


D 


E 


Parallel 


0.5 


-0.15 


0.3 


-0.02223 


-0.02223 


Perpendicular 


0.5 


-0.25 


0.5 


-0.04167 


0.04167 



Table 1 . Fit parameters A, B, C, D, and E in Eq. (4) in 
the case of the stripes oriented parallel to the force, 
and perpendicular to the driving force. 



analytically by Bazant and Vinogradova ( 12), to be 

b = b\\ cos 2 + b ± sin 2 . (5) 

We find that the slip length obtained in our simulations 
using the tilted force agrees with the theoretical prediction 
of Eq. §5§ with a relative error of only 10 _3 %. In Fig. [2] 
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slip parameter 

Figure 2. Slip length of a striped wall as a function 
of the slip parameter of the stripes for four different 
orientations of the stripes with respect to the direction 
of the driving force: parallel (<d = 0°), perpendicular 
(tf = 90°), diagonal (tf = 45°), and tilted by 60 degrees. 
The condition Ci + C2 = 1 is maintained, leading to 
symmetry about C = 0.5. 



we show the dependence of the effective slip length on the 
orientation of the force and on the choice of slip parameters 
applied on the stripes. For studying the angle dependence, 
the stripes have been constrained by the condition Ci + C2 = 
1. For (1 = (2 = 0.5, the case of a homogeneous wall is 
restored and therefore, the slip length is independent on the 
angle $ between the force and the direction of the extension 
of the stripes. If £i and (2 are chosen differently, the effective 
slip depends on the orientation. For $ = 0°, i.e. if the force 
is aligned parallel to the stripes, a maximal slip length is 
obtained. If the stripes are oriented perpendicular to the 
direction of the force, a minimal slip length is obtained. Out 
of those extreme cases the slip length for the orientations in 
between can be calculated from Eq. ([5]). 
We also find that in simulations with the force acting in 
diagonal direction, the fluid tends to follow the stripes on 
the wall. The velocity of the fluid shows a small component 
which is perpendicular to the force. This component can be 
understood as a consequence of the tensorial nature of the 
slip as proposed by Bazant and Vinogradova (fT2l). In our 
simulations, the transverse flux is four orders of magnitude 
less than the velocity obtained in the direction of the force. 
Even though the effect is small, these simulations confirm 
that walls with stripes of varying slip length can be used in 
mixing devices to cause vortices in micro and nano flows. 

3.2. Continuously varying Striped Walls 

In this section, the slip parameter ( is made to vary as a 
continuous periodic function along the wall nodes in the y 
and z directions, in the following manner: 

C = C + C-cos(k-x). (6) 



where the wave vector k defines the frequency and direction 
of the variation of slip parameter (. We use an amplitude 
( = 0.5 and a mean slip parameter of ( = 0.5. 




Figure 3. Walls are designed having stripes with con- 
tinuously varying slip parameter C, as per Eq. (6). The 
colors correspond to the velocity projected in the xy- 
plane scaled by a factor of 3000. 

We run simulations having periodicity of the slip parameter 
as A = ^r>/2 lattice nodes, giving rise to stripes as shown 
in Fig. [3] The wave vector points in the diagonal direction 
on each of the boundary planes and the magnitude is chosen 
such that along each side two stripes fit between the edge 
nodes. 




Figure 4. The velocity along the x boundary in the case 
of variation in the slip parameter along the wall nodes. 
The transverse velocity component is color coded after 
scaling by a factor of 3000. 

The influence of such a wall design can be understood by 
plotting the velocity vectors along the wall, as shown in 
Fig. [4] It is clearly seen that as ( tends to 1, the velocity 
along such nodes is increased, and as ( tends to 0, the 
velocity along such nodes is reduced. As a result the fluid 
arriving in a region where the slip length decreases has 
to go aside. Hence the occurrence of a vortice in the 
plane perpendicular to the direction of the force driven 



flow may be observed. Due to the varying slip conditions a 
varying direction of the velocity occurs along the z direction. 
However, if averaged over all lattice planes a net flux in a 
vortex remains, as shown in Fig. [5] 

In order to better quantify the influence of the pattern on the 
velocity profile, we calculate the z-component of the curl of 
the velocity field. As shown in Fig.[6j close to the surface, 
there are distortions by the pattern, but in the center of the 
channel a homogeneous rotation of the velocity field about 
the axis of the channel is observed. 

If the z-component of the curl is averaged along the z- 
direction, indeed a net rotation which is constant in large 
areas of the cross section is identified. In Fig.|7]this is plotted 
for both cases of small and large periodicity of the pattern. 
One can see that with the small period of the pattern the 
rotation is more homogeneous, whereas with a large period 
the distortion of the flow is stronger, and the curl obtained 
in the center is slightly enhanced. 




Figure 5. An averaged projection of velocity vectors of 
the xy-p\ane in the case of variation in the slip param- 
eter along the wall nodes, along the entire z-direction 
domain. It shows the existence of a homogeneous 
rotation about the centre. The velocity at the centre is 
aligned with the channel, whereas everywhere else a 
net circular is seen. 

Our results confirm previous understanding that at a de- 
creased angle (stripes along the driving force) the slip 
length increases. We further understand that if the stripes 
are oriented parallel to the force, the fluid can flow along 
the stripes and form lamellae of fluid with different slip 
according to the local surface properties on the respective 
stripe. 

4. Conclusion 

In this paper we have applied the boundary condition of 
Ref. (7) to patterned surfaces. We have studied the depen- 
dence of the effective slip length in channels with textured 
surfaces of stripes of equal width, with alternating local slip 
length. We find that the effective slip length can be expressed 
as function of the individual slip lengths on the stripes. 




Figure 6. A 3D plot of the z-component of the curl of 
the velocity field in the case of variation in £ along the 
wall nodes, shows the distortions of the flow field by 
the pattern on the wall and reflects the homogeneous 
rotation about the centre again. The value calculated 
for the z-component of the curl is scaled by a factor of 
10 4 for plotting. 




Figure 7. The z-component of the curl of the velocity 
field averaged along the z-direction in the case of vari- 
ation in C along the wall nodes. Here the homogeneous 
rotation at the center of the channel is reflected by the 
constant value of the curl. The value calculated for the 
z-component of the curl is scaled by a factor of 10 4 for 
plotting. 



The tensorial nature of the slip can be exploited to make 
the fluid flow follow a given direction. We have tested 
this in simulations where the flow through a rectangular 
channel patterned with diagonal stripes is simulated. The 
flow follows the stripes and so a vortex is induced in the flow 
field. This can be used to construct a micro mixer device. 
Our LB simulations can help to optimize the patterns to put 
on the surfaces and to study the dependence of the flow field 
on various parameters, such as the period of the pattern or 
the wetting contrast between two stripes. 
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